Setup chunk
Setup reticulate
knitr::opts_chunk$set(fig.width = 8)
knitr::opts_knit$set(root.dir = normalizePath(".."))
knitr::opts_knit$get("root.dir")
[1] "/nas/groups/treutlein/USERS/tomasgomes/projects/liver_regen"
Load libraries
library(reticulate)
knitr::knit_engines$set(python = reticulate::eng_python)
py_available(initialize = FALSE)
[1] FALSE
use_python(Sys.which("python"))
py_config()
python: /home/tpires/bin/miniconda3/bin/python
libpython: /home/tpires/bin/miniconda3/lib/libpython3.8.so
pythonhome: /home/tpires/bin/miniconda3:/home/tpires/bin/miniconda3
version: 3.8.3 (default, May 19 2020, 18:47:26) [GCC 7.3.0]
numpy: /home/tpires/bin/miniconda3/lib/python3.8/site-packages/numpy
numpy_version: 1.22.1
NOTE: Python version was forced by RETICULATE_PYTHON
Other functions
library(Seurat)
Attaching SeuratObject
library(ggplot2)
library(plyr)
library(dplyr)
Attaching package: ‘dplyr’
The following objects are masked from ‘package:plyr’:
arrange, count, desc, failwith, id, mutate, rename, summarise, summarize
The following objects are masked from ‘package:stats’:
filter, lag
The following objects are masked from ‘package:base’:
intersect, setdiff, setequal, union
library(RColorBrewer)
Load data
# downsample taking the minimum median of UMI counts as a reference
dsSeurat = function(srat, variable, assay = "SCT", numis = "nCount_SCT", seed = 1){
# input Seurat object and variable to do the downsample based on
# return Seurat object
set.seed(seed)
# downsample taking the minimum median of UMI counts as a reference
fact = srat@meta.data[,numis]/min(tapply(srat@meta.data[,numis], srat@meta.data[,variable], median))
fact = ifelse(fact<=1, 0, log2(fact))
thinned_counts = seqgendiff::thin_lib(as.matrix(srat@assays[[assay]]@counts),
thinlog2 = fact, relative = T)
rownames(thinned_counts$mat) = rownames(srat@assays[[assay]]@data)
colnames(thinned_counts$mat) = colnames(srat@assays[[assay]]@data)
srat@assays$SCT@counts = Matrix::Matrix(thinned_counts$mat, sparse = T)
srat@assays$SCT@data = log1p(srat@assays[[assay]]@counts)
return(srat)
}
Get downsampled hepatocytes
f = "results/cirrhosis/sub_allcells_css.RDS"
if(!file.exists(f)){
allcells_css = readRDS(file = "data/processed/allcells_css.RDS")
end_cells = readRDS(file = "results/endothelial/only_end_cells_zon.RDS")
hep_cells = readRDS(file = "results/zonation_cond/hep_cells_zonation_rank.RDS")
imm_cells = readRDS(file = "results/immune/all_imm_cells.RDS")
endpop_df = data.frame(row.names = colnames(end_cells),
"subpops" = end_cells@meta.data$endo_simp)
immpop_df = data.frame(row.names = colnames(imm_cells),
"subpops" = imm_cells@meta.data$immune_annot)
heppop_df = lapply(hep_cells, function(x) cbind(colnames(x), as.character(x$zonation_int)))
heppop_df = Reduce(rbind, heppop_df)
heppop_df = data.frame(row.names = heppop_df[,1],
subpops = factor(heppop_df[,2]))
levels(heppop_df$subpops) = c("(-0.00099,0.333]" = "Hepatocytes_Z1",
"(0.333,0.667]" = "Hepatocytes_Z2",
"(0.667,1]" = "Hepatocytes_Z3")
heppop_df$subpops = as.character(heppop_df$subpops)
subpop_df = rbind(endpop_df, immpop_df, heppop_df)
subpop_df$subpops[subpop_df$subpops=="Cycling cells"] = "Dividing endothelial cells"
# add subpop metadata
allcells_css = AddMetaData(allcells_css, subpop_df)
allcells_css$subpops[is.na(allcells_css$subpops)] = allcells_css$allcells_major[is.na(allcells_css$subpops)]
allcells_css$subpops[allcells_css$subpops=="Hepatocyte-Monocyte interaction"] = "Doublets"
# the cells in this object named "Hepatocytes", "Endothelial cells", and "Doublets" have to be removed
## the first two are cells that didn't pass the hep and end analysis
sub_allcells_css = allcells_css[,!(allcells_css$subpops %in% c("Hepatocytes", "Endothelial cells",
"Doublets"))]
allcells_css$subpops[allcells_css$subpops %in% c("Hepatocytes", "Endothelial cells","Doublets")] = "Not annotated"
sub_allcells_css$anno_cond = paste0(sub_allcells_css$subpops, "_", sub_allcells_css$Condition)
saveRDS(sub_allcells_css, file = f)
} else{
sub_allcells_css = readRDS(f)
}
Get cell type markers
f = "results/cirrhosis/ds_hep_srat.RDS"
if(!file.exists(f)){
ds_hep_srat = dsSeurat(sub_allcells_css[,grepl("Hepatoc", sub_allcells_css$subpops)], "Condition")
saveRDS(ds_hep_srat, file = f)
} else{
ds_hep_srat = readRDS(f)
}
Load Ramachandran data
cell_type_mk = readRDS(file = "results/cond_effect_v2/cell_type_mk_more.RDS")
Load Camp data
f = "results/cirrhosis/cirr_data_renorm.RDS"
if(!file.exists(f)){
load("../../data/published/Ramachandran_liver/tissue.rdata")
cirr_data = CreateSeuratObject(counts = tissue@raw.data, meta.data = tissue@meta.data)
# Renormalise, since the original data doesn't look ok
cirr_data = suppressWarnings(SCTransform(cirr_data, do.correct.umi = T, verbose = F,
vars.to.regress=c("dataset", "nCount_RNA"),
variable.features.rv.th = 1, seed.use = 1,
return.only.var.genes = F, variable.features.n = NULL))
cirr_data$anno_cond = paste0(cirr_data$annotation_indepth, "_", cirr_data$condition)
saveRDS(cirr_data, file = f)
} else{
cirr_data = readRDS(f)
}
Get DE genes between conditions
f = "results/cirrhosis/mk_ct_CvsH.RDS"
if(!file.exists(f)){
cirr_data = SetIdent(cirr_data, value = "annotation_indepth")
mk_ct_list = list()
condct = table(cirr_data$annotation_indepth, cirr_data$condition)
condct = condct[apply(condct, 1, function(x) all(x>=5)),]
for(ct in rownames(condct)){
mk_ct_list[[ct]] = FindMarkers(cirr_data, ident.1 = "Cirrhotic", ident.2 = "Uninjured",
assay = "SCT", group.by = "condition", pseudocount.use = 0.1,
subset.ident = ct)
}
saveRDS(mk_ct_list, file = f)
} else{
mk_ct_list = readRDS(f)
}
f = "results/cirrhosis/mk_li_CvsH.RDS"
if(!file.exists(f)){
cirr_data = SetIdent(cirr_data, value = "annotation_lineage")
mk_li_list = list()
condct = table(cirr_data$annotation_lineage, cirr_data$condition)
condct = condct[apply(condct, 1, function(x) all(x>=5)),]
for(ct in rownames(condct)){
mk_li_list[[ct]] = FindMarkers(cirr_data, ident.1 = "Cirrhotic", ident.2 = "Uninjured",
assay = "SCT", group.by = "condition",
pseudocount.use = 0.1, subset.ident = ct)
}
saveRDS(mk_li_list, file = f)
} else{
mk_li_list = readRDS(f)
}
f = "results/cirrhosis/mk_cond.RDS"
if(!file.exists(f)){
mk_cond = FindMarkers(cirr_data, ident.1 = "Cirrhotic", ident.2 = "Uninjured", assay = "SCT",
group.by = "condition", pseudocount.use = 0.1)
saveRDS(mk_cond, file = f)
} else{
mk_cond = readRDS(f)
}
Markers for each MPs cell population
f = "results/cirrhosis/mk_ct_CvsH.RDS"
if(!file.exists(f)){
cirr_data = SetIdent(cirr_data, value = "annotation_indepth")
mk_ct_list = list()
condct = table(cirr_data$annotation_indepth, cirr_data$condition)
condct = condct[apply(condct, 1, function(x) all(x>=5)),]
for(ct in rownames(condct)){
mk_ct_list[[ct]] = FindMarkers(cirr_data, ident.1 = "Cirrhotic", ident.2 = "Uninjured",
assay = "SCT", group.by = "condition", pseudocount.use = 0.1,
subset.ident = ct)
}
saveRDS(mk_ct_list, file = f)
} else{
mk_ct_list = readRDS(f)
}
f = "results/cirrhosis/mk_li_CvsH.RDS"
if(!file.exists(f)){
cirr_data = SetIdent(cirr_data, value = "annotation_lineage")
mk_li_list = list()
condct = table(cirr_data$annotation_lineage, cirr_data$condition)
condct = condct[apply(condct, 1, function(x) all(x>=5)),]
for(ct in rownames(condct)){
mk_li_list[[ct]] = FindMarkers(cirr_data, ident.1 = "Cirrhotic", ident.2 = "Uninjured",
assay = "SCT", group.by = "condition",
pseudocount.use = 0.1, subset.ident = ct)
}
saveRDS(mk_li_list, file = f)
} else{
mk_li_list = readRDS(f)
}
f = "results/cirrhosis/mk_cond.RDS"
if(!file.exists(f)){
mk_cond = FindMarkers(cirr_data, ident.1 = "Cirrhotic", ident.2 = "Uninjured", assay = "SCT",
group.by = "condition", pseudocount.use = 0.1)
saveRDS(mk_cond, file = f)
} else{
mk_cond = readRDS(f)
}
Markers for each cell population (overall and within compartment)
f = "results/cirrhosis/mk_cirr_mps.RDS"
if(!file.exists(f)){
cirr_data = SetIdent(cirr_data, value = "annotation_indepth")
mk_cirr_mps = FindAllMarkers(cirr_data[,grepl("MPs", cirr_data$annotation_indepth)],
assay = "SCT", pseudocount.use = 0.1, only.pos = T,
logfc.threshold = 0.15)
mk_cirr_mps = mk_cirr_mps[mk_cirr_mps$p_val_adj<=0.05,]
mk_cirr_mps = mk_cirr_mps[order(mk_cirr_mps$avg_log2FC, decreasing = T),]
saveRDS(mk_cirr_mps, file = f)
} else{
mk_cirr_mps = readRDS(f)
}
Scoring global signature
f = "results/cirrhosis/mk_ct_all.RDS"
if(!file.exists(f)){
mk_ct_all = list()
for(comp in c("annotation_lineage", "annotation_indepth")){
cirr_data = SetIdent(cirr_data, value = comp)
mk_ct = FindAllMarkers(cirr_data, assay = "SCT", pseudocount.use = 0.1, only.pos = T,
logfc.threshold = 0.15)
mk_ct = mk_ct[mk_ct$p_val_adj<=0.05,]
mk_ct = mk_ct[order(mk_ct$avg_log2FC, decreasing = T),]
mk_ct_all[[comp]] = mk_ct
if(comp=="annotation_indepth"){
l_sub = list()
for(lin in unique(cirr_data$annotation_lineage)){
print(lin)
nct = sum(table(cirr_data$annotation_indepth[cirr_data$annotation_lineage==lin])>0)
if(nct>1){
mk_sub = FindAllMarkers(cirr_data[,cirr_data$annotation_lineage==lin], only.pos = T,
assay = "SCT", pseudocount.use = 0.1, logfc.threshold = 0.3)
mk_sub = mk_sub[mk_sub$p_val_adj<=0.05,]
mk_sub = mk_sub[order(mk_sub$avg_log2FC, decreasing = T),]
l_sub[[lin]] = mk_sub
}
}
mk_ct_all[["perlineage"]] = Reduce(rbind, l_sub)
}
}
saveRDS(mk_ct_all, file = f)
} else{
mk_ct_all = readRDS(f)
}
Scoring by signature determined within lineage (downsamples data to match conditions)
cir_sig_all = list("cirrhotic_all" = rownames(mk_cond)[mk_cond$avg_log2FC>=0.3 & mk_cond$pct.1>0.2 &
mk_cond$pct.1-mk_cond$pct.2>0.1])
sub_allcells_css = AddModuleScore(sub_allcells_css, features = cir_sig_all, nbin = 20, seed = 1,
ctrl = 500, name = "cirrhotic_all")
for(lin in unique(sub_allcells_css$allcells_major)){
ids = unique(sub_allcells_css@active.ident[sub_allcells_css$allcells_major==lin])
plt = VlnPlot(sub_allcells_css, group.by = "subpops", split.by = "Condition",
features = "cirrhotic_all1", idents = ids, ncol = 1)+theme(legend.position = "bottom")
print(plt)
}
The default behaviour of split.by has changed.
Separate violin plots are now plotted side-by-side.
To restore the old behaviour of a single split violin,
set split.plot = TRUE.
This message will be shown once per session.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
Warning: Groups with fewer than two data points have been dropped.
List, for each cell type, which signatures could make sense
feats = colnames(sub_allcells_dshep@meta.data)[col_sigs]
sig_ass_list = lapply(unique(sub_allcells_dshep@meta.data$subpops), function(x) c())
names(sig_ass_list) = unique(sub_allcells_dshep@meta.data$subpops)
for(ct in names(sig_ass_list)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], "Cirr_li_all")
if(grepl("NK ", ct) | grepl("Treg", ct) | grepl("ILC", ct) |
grepl("T cell", ct) | grepl("MAIT", ct) | grepl("TRM", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("ILCs", feats) | grepl("Tcells", feats)])
} else if(grepl("DC", ct) | grepl("Kupffer", ct) | grepl("Monocyte", ct) | grepl("Macrophage", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("MPs", feats) | grepl("DC", feats)])
} else if(grepl("Hepato", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], "Cirr_li_Hepatocytes")
} else if(grepl("Cholangiocytes", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Cholangiocytes", feats) ])
} else if(grepl("B cell", ct) | grepl("Plasma", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Bcells", feats)])
} else if(grepl("Stellate", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Mesenchyme", feats) | grepl("Myofibroblast", feats)])
} else if(grepl("EC", ct) | grepl("endothelial", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Endothelia", feats)])
}
}
Testing for significance and individual plotting
feats = colnames(sub_allcells_dshep@meta.data)[col_sigs]
sig_ass_list = lapply(unique(sub_allcells_dshep@meta.data$subpops), function(x) c())
names(sig_ass_list) = unique(sub_allcells_dshep@meta.data$subpops)
for(ct in names(sig_ass_list)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], "Cirr_li_all")
if(grepl("NK ", ct) | grepl("Treg", ct) | grepl("ILC", ct) |
grepl("T cell", ct) | grepl("MAIT", ct) | grepl("TRM", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("ILCs", feats) | grepl("Tcells", feats)])
} else if(grepl("DC", ct) | grepl("Kupffer", ct) | grepl("Monocyte", ct) | grepl("Macrophage", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("MPs", feats) | grepl("DC", feats)])
} else if(grepl("Hepato", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], "Cirr_li_Hepatocytes")
} else if(grepl("Cholangiocytes", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Cholangiocytes", feats) ])
} else if(grepl("B cell", ct) | grepl("Plasma", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Bcells", feats)])
} else if(grepl("Stellate", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Mesenchyme", feats) | grepl("Myofibroblast", feats)])
} else if(grepl("EC", ct) | grepl("endothelial", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Endothelia", feats)])
}
}
Check top hits
pvaldf_filt = pvaldf[pvaldf$fdr<=0.05,]
sigfc_plts = list()
for(cc in c("healthy_v_embolised", "healthy_v_regenerating")){
subdf = pvaldf_filt[pvaldf_filt$L3==cc,]
sigfc_plts[[cc]] = list()
for(ct in unique(pvaldf_filt$L1)){
subsubdf = subdf[subdf$L1==ct,]
subsubdf = subsubdf[order(subsubdf$fc, decreasing = F),]
subsubdf$L2 = factor(subsubdf$L2, levels = unique(subsubdf$L2))
sigfc_plts[[cc]][[ct]] = ggplot(subsubdf, aes(x = L2, y = fc))+
geom_col()+
geom_hline(yintercept = c(-.2, .2), linetype = "dashed")+
theme_classic()+
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1))
}
}
Calculate signatures for each cell type
sub_allcells_dshep@meta.data = sub_allcells_dshep@meta.data[,1:118]
# get top markers for cell types
cttop = list()
topn = 50
for(i in unique(mk_ct_all$perlineage$cluster)){
subdf = mk_ct_all$perlineage[mk_ct_all$perlineage$cluster==i & mk_ct_all$perlineage$p_val_adj<=0.05,]
cttop[[i]] = subdf[order(subdf$avg_log2FC, decreasing = T),"gene"][1:topn]
cttop = lapply(cttop, function(x) x[!is.na(x)])
}
# get signatures
sub_allcells_dshep = AddModuleScore(sub_allcells_dshep, features = cttop, nbin = 15,
seed = 1, ctrl = 500, name = "ct_sig_")
col_sigs = grepl("ct_sig_", colnames(sub_allcells_dshep@meta.data))
colnames(sub_allcells_dshep@meta.data)[col_sigs] = paste0("ct_sig_", names(cttop))
List, for each cell type, which signatures could make sense
feats = colnames(sub_allcells_dshep@meta.data)[startsWith(colnames(sub_allcells_dshep@meta.data), "ct_sig_")]
sig_ass_list = lapply(unique(sub_allcells_dshep@meta.data$subpops), function(x) c())
names(sig_ass_list) = unique(sub_allcells_dshep@meta.data$subpops)
for(ct in names(sig_ass_list)){
if(grepl("NK ", ct) | grepl("Treg", ct) | grepl("ILC", ct) |
grepl("T cell", ct) | grepl("MAIT", ct) | grepl("TRM", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("ILCs", feats) | grepl("Tcells", feats)])
} else if(grepl("DC", ct) | grepl("Kupffer", ct) | grepl("Monocyte", ct) | grepl("Macrophage", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("MPs", feats) | grepl("DC", feats)])
} else if(grepl("Hepato", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], "Cirr_li_Hepatocytes")
} else if(grepl("Cholangiocytes", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Cholangiocytes", feats) ])
} else if(grepl("B cell", ct) | grepl("Plasma", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Bcells", feats)])
} else if(grepl("Stellate", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Mesenchyme", feats) | grepl("Myofibroblast", feats)])
} else if(grepl("EC", ct) | grepl("endothelial", ct)){
sig_ass_list[[ct]] = c(sig_ass_list[[ct]], feats[grepl("Endothelia", feats)])
}
}
Plot signatures for groups of cell types
ct_vars = grepl("ct_sig_", colnames(sub_allcells_dshep@meta.data))
plot_df = sub_allcells_dshep@meta.data[,c("subpops", "Condition",
colnames(sub_allcells_dshep@meta.data)[ct_vars])]
plot_df = reshape2::melt(plot_df)
plot_ct_sigs = list()
for(gr in c("Kupffer cells", "CD8 ab-T cells 1", "EC non-LSEC", "Cholangiocytes")){
ct_use = names(sig_ass_list)[unlist(lapply(sig_ass_list,
function(x) length(intersect(x,sig_ass_list[[gr]]))==length(x)))]
sub_plot_df = plot_df[plot_df$subpops %in% ct_use & plot_df$variable %in% sig_ass_list[[gr]],]
sub_plot_df$variable = substr(sub_plot_df$variable , 8, 200)
mean_df = reshape2::melt(tapply(sub_plot_df$value, sub_plot_df$variable, mean))
colnames(mean_df)[1] = "variable"
plot_ct_sigs[[gr]] = ggplot(sub_plot_df, aes(x = value, y = subpops, fill = subpops))+
facet_wrap(~variable, ncol = length(mean_df$variable), scales = "free_x")+
geom_violin(size = 0.2)+
geom_vline(mapping = aes(xintercept = value), data = mean_df, linetype = "dashed", size = 0.4)+
labs(x = "Signature score", y = "Subpopulation")+
theme_classic()+
theme(legend.position = "none",
axis.title = element_text(size = 8),
axis.text.y = element_text(size = 7, colour = "black"),
strip.text = element_text(size = 7, colour = "black"),
axis.text.x = element_text(size = 6.6, colour = "black"))
}
Testing for significance and individual plotting
cond_mat = combn(unique(sub_allcells_dshep$Condition), 2)
pval_list = list()
diff_list = list()
plots_list = list()
for(s in names(sig_ass_list)){
pval_list[[s]] = list()
diff_list[[s]] = list()
plots_list[[s]] = list()
for(f in sig_ass_list[[s]]){
pval_list[[s]][[f]] = list()
diff_list[[s]][[f]] = list()
plots_list[[s]][[f]] = VlnPlot(sub_allcells_dshep, group.by = "subpops", idents = s, features = f,
split.by = "Condition", ncol = 1)+theme(legend.position = "bottom")
# rescaling (not for plotting only testing)
if(!(paste0(f, "_resc") %in% colnames(sub_allcells_dshep@meta.data))){
sub_allcells_dshep@meta.data[,paste0(f, "_resc")] = scales::rescale(sub_allcells_dshep@meta.data[,f],
to = c(0, max(scales::rescale(sub_allcells_dshep@meta.data[,f]))))
}
for(i in 1:ncol(cond_mat)){
test = wilcox.test(sub_allcells_dshep@meta.data[sub_allcells_dshep$subpops==s &
sub_allcells_dshep$Condition==cond_mat[1,i],paste0(f, "_resc")],
sub_allcells_dshep@meta.data[sub_allcells_dshep$subpops==s &
sub_allcells_dshep$Condition==cond_mat[2,i],paste0(f, "_resc")],
exact = F)
id = paste0(cond_mat[1,i], "_v_", cond_mat[2,i])
pval_list[[s]][[f]][[id]] = test$p.value
d = mean(sub_allcells_dshep@meta.data[sub_allcells_dshep$subpops==s &
sub_allcells_dshep$Condition==cond_mat[1,i],paste0(f, "_resc")])/mean(sub_allcells_dshep@meta.data[sub_allcells_dshep$subpops==s &
sub_allcells_dshep$Condition==cond_mat[2,i],paste0(f, "_resc")])
diff_list[[s]][[f]][[id]] = log2(d)
}
}
}
pvaldf_ct = reshape2::melt(pval_list)[,4:1]
pvaldf_ct$fdr = p.adjust(pvaldf_ct$value, method = "fdr")
diffdf_ct = reshape2::melt(diff_list)[,4:1]
pvaldf_ct$fc = diffdf_ct$value
Calculate cirrhotic signatures
sigs_endo_cond = cir_sigsall[c("Endothelia (1)","Endothelia (2)","Endothelia (3)",
"Endothelia (5)","Endothelia (6)", "Endothelia (7)")]
end_cells@meta.data = end_cells@meta.data[,!grepl("cirr_cond_", colnames(end_cells@meta.data))]
end_cells = AddModuleScore(end_cells, features = sigs_endo_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_cond_")
col_sigs = grepl("cirr_cond_", colnames(end_cells@meta.data))
colnames(end_cells@meta.data)[col_sigs] = paste0("cirr_cond_", names(sigs_endo_cond))
colnames(end_cells@meta.data) = gsub(" (", "_", colnames(end_cells@meta.data), fixed = T)
colnames(end_cells@meta.data) = gsub(")", "", colnames(end_cells@meta.data), fixed = T)
sigs_mono_cond = cir_sigsall[c("MPs (1)","MPs (2)","MPs (3)","MPs (4)","MPs (5)","MPs (6)",
"MPs (7)","MPs (8)","MPs (9)",
"Cycling MPs (1)","Cycling MPs (3)")]
mon_cells@meta.data = mon_cells@meta.data[,!grepl("cirr_cond_", colnames(mon_cells@meta.data))]
mon_cells = AddModuleScore(mon_cells, features = sigs_mono_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_cond_")
col_sigs = grepl("cirr_cond_", colnames(mon_cells@meta.data))
colnames(mon_cells@meta.data)[col_sigs] = paste0("cirr_cond_", names(sigs_mono_cond))
colnames(mon_cells@meta.data) = gsub(" (", "_", colnames(mon_cells@meta.data), fixed = T)
colnames(mon_cells@meta.data) = gsub(")", "", colnames(mon_cells@meta.data), fixed = T)
sigs_hepa_cond = cir_sigsall[c("Hepatocytes")]
ds_hep_srat@meta.data = ds_hep_srat@meta.data[,!grepl("cirr_cond_",
colnames(ds_hep_srat@meta.data))]
ds_hep_srat = AddModuleScore(ds_hep_srat, features = sigs_hepa_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_cond_")
col_sigs = grepl("cirr_cond_", colnames(ds_hep_srat@meta.data))
colnames(ds_hep_srat@meta.data)[col_sigs] = paste0("cirr_cond_", names(sigs_hepa_cond))
Calculate cell type signatures
sigs_endo_cond = cir_sigsall[c("Endothelia (1)","Endothelia (2)","Endothelia (3)",
"Endothelia (5)","Endothelia (6)", "Endothelia (7)")]
end_cells@meta.data = end_cells@meta.data[,!grepl("cirr_cond_", colnames(end_cells@meta.data))]
end_cells = AddModuleScore(end_cells, features = sigs_endo_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_cond_")
col_sigs = grepl("cirr_cond_", colnames(end_cells@meta.data))
colnames(end_cells@meta.data)[col_sigs] = paste0("cirr_cond_", names(sigs_endo_cond))
colnames(end_cells@meta.data) = gsub(" (", "_", colnames(end_cells@meta.data), fixed = T)
colnames(end_cells@meta.data) = gsub(")", "", colnames(end_cells@meta.data), fixed = T)
sigs_mono_cond = cir_sigsall[c("MPs (1)","MPs (2)","MPs (3)","MPs (4)","MPs (5)","MPs (6)",
"MPs (7)","MPs (8)","MPs (9)",
"Cycling MPs (1)","Cycling MPs (3)")]
mon_cells@meta.data = mon_cells@meta.data[,!grepl("cirr_cond_", colnames(mon_cells@meta.data))]
mon_cells = AddModuleScore(mon_cells, features = sigs_mono_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_cond_")
col_sigs = grepl("cirr_cond_", colnames(mon_cells@meta.data))
colnames(mon_cells@meta.data)[col_sigs] = paste0("cirr_cond_", names(sigs_mono_cond))
colnames(mon_cells@meta.data) = gsub(" (", "_", colnames(mon_cells@meta.data), fixed = T)
colnames(mon_cells@meta.data) = gsub(")", "", colnames(mon_cells@meta.data), fixed = T)
sigs_hepa_cond = cir_sigsall[c("Hepatocytes")]
ds_hep_srat@meta.data = ds_hep_srat@meta.data[,!grepl("cirr_cond_",
colnames(ds_hep_srat@meta.data))]
ds_hep_srat = AddModuleScore(ds_hep_srat, features = sigs_hepa_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_cond_")
col_sigs = grepl("cirr_cond_", colnames(ds_hep_srat@meta.data))
colnames(ds_hep_srat@meta.data)[col_sigs] = paste0("cirr_cond_", names(sigs_hepa_cond))
Plot signatures - endothelial
sigs_endo_cond = cttop[c("Endothelia (1)","Endothelia (2)","Endothelia (3)","Endothelia (4)",
"Endothelia (5)","Endothelia (6)", "Endothelia (7)")]
end_cells@meta.data = end_cells@meta.data[,!grepl("cirr_ct_", colnames(end_cells@meta.data))]
end_cells = AddModuleScore(end_cells, features = sigs_endo_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_ct_")
col_sigs = grepl("cirr_ct_", colnames(end_cells@meta.data))
colnames(end_cells@meta.data)[col_sigs] = paste0("cirr_ct_", names(sigs_endo_cond))
colnames(end_cells@meta.data) = gsub(" (", "_", colnames(end_cells@meta.data), fixed = T)
colnames(end_cells@meta.data) = gsub(")", "", colnames(end_cells@meta.data), fixed = T)
sigs_mono_cond = cttop[c("MPs (1)","MPs (2)","MPs (3)","MPs (4)","MPs (5)","MPs (6)",
"MPs (7)","MPs (8)","MPs (9)", "Cycling MPs (1)",
"Cycling MPs (2)","Cycling MPs (3)","Cycling MPs (4)")]
mon_cells@meta.data = mon_cells@meta.data[,!grepl("cirr_ct_", colnames(mon_cells@meta.data))]
mon_cells = AddModuleScore(mon_cells, features = sigs_mono_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_ct_")
Warning: The following features are not present in the object: PRRG3, ATXN8OS, not searching for symbol synonyms
col_sigs = grepl("cirr_ct_", colnames(mon_cells@meta.data))
colnames(mon_cells@meta.data)[col_sigs] = paste0("cirr_ct_", names(sigs_mono_cond))
colnames(mon_cells@meta.data) = gsub(" (", "_", colnames(mon_cells@meta.data), fixed = T)
colnames(mon_cells@meta.data) = gsub(")", "", colnames(mon_cells@meta.data), fixed = T)
sigs_hepa_cond = cttop[c("Hepatocytes")]
ds_hep_srat@meta.data = ds_hep_srat@meta.data[,!grepl("cirr_ct_",
colnames(ds_hep_srat@meta.data))]
ds_hep_srat = AddModuleScore(ds_hep_srat, features = sigs_hepa_cond, nbin = 30,
seed = 1, ctrl = 500, name = "cirr_ct_")
col_sigs = grepl("cirr_ct_", colnames(ds_hep_srat@meta.data))
colnames(ds_hep_srat@meta.data)[col_sigs] = paste0("cirr_ct_", names(sigs_hepa_cond))
Plot signatures - monocytes
VlnPlot(end_cells, group.by = "endo_simp",
features = c("cirr_cond_Endothelia_1", "cirr_cond_Endothelia_2",
"cirr_cond_Endothelia_3", "cirr_cond_Endothelia_5",
"cirr_cond_Endothelia_6", "cirr_cond_Endothelia_7"),
split.by = "Condition", ncol = 1)
VlnPlot(end_cells, group.by = "endo_simp",
features = c("cirr_ct_Endothelia_1", "cirr_ct_Endothelia_2", "cirr_ct_Endothelia_3",
"cirr_ct_Endothelia_4", "cirr_ct_Endothelia_5", "cirr_ct_Endothelia_6",
"cirr_ct_Endothelia_7"),
split.by = "Condition", ncol = 1)
Plot signatures - hepatocytes
VlnPlot(mon_cells, group.by = "mono_annot",
features = c("cirr_cond_MPs_1", "cirr_cond_MPs_2", "cirr_cond_MPs_3", "cirr_cond_MPs_4",
"cirr_cond_MPs_5", "cirr_cond_MPs_6","cirr_cond_MPs_7", "cirr_cond_MPs_8",
"cirr_cond_MPs_9", "cirr_cond_Cycling MPs_1", "cirr_cond_Cycling MPs_3"),
split.by = "Condition", ncol = 1)